#include <stdio.h>
#include <gmp.h>
#include <math.h>
#include <pthread.h>

#define IT 1000000000
#define PREC 1000

// dados
typedef struct {
	mpf_t a_prev, a_next;
	mpf_t b_prev, b_next;
	mpf_t p_prev, p_next;
}data;

void *calc_a (void *param);
void *calc_b (void *param);
void *calc_p (void *param);


int main(int argc, char *argv[]) {

	//Verefica os dados de entrada 
	if (argc != 2) {
		printf("usage: ./borwein (0 = sequential / 1 = parallel)\n");
		return -1;
	}

	// atribui PREC como precisão
	mpf_set_default_prec(PREC);

	// declara dados
	data values;

	// inicializa variaveis
	mpf_init(values.a_prev); mpf_init(values.a_next);
	mpf_init(values.b_prev); mpf_init(values.b_next);
	mpf_init(values.p_prev); mpf_init(values.p_next);

	// Variavel temporaria para calculos
	mpf_t tmp;
	mpf_init(tmp);

	double d;
	int i;

	// Atribui valores iniciais para a_next e p_next

	// a_next = sqrt(2)
	d = sqrt(2);
	mpf_set_d(values.a_next, d);
	
	// p_next = 2 + sqrt(2)
	d = 2 + sqrt(2);
	mpf_set_d(values.p_next, d);

	// INICIO DO PROCESSO EM SEQUENCIAL
	if(argv[1][0] == '0') {


		// Inicia interação
		for(i=0; i<IT; ++i) {

			// CALCULA A_NEXT
			// a_prev = a_next
			mpf_set(values.a_prev, values.a_next);
			// tmp = sqrt(a_prev)
			mpf_sqrt(tmp, values.a_prev);
			// a_next = 1 / tmp
			mpf_ui_div(values.a_next, 1, tmp);
			// a_next = a_next + tmp
			mpf_add(values.a_next, values.a_next, tmp);
			// a_next = a_next / 2
			mpf_div_ui(values.a_next, values.a_next, 2);

			// CALCULA B_NEXT
			// b_prev = b_next
			mpf_set(values.b_prev, values.b_next);
			// b_next = b_prev + 1
			mpf_add_ui(values.b_next, values.b_prev, 1);
			// b_next = b_next * tmp
			mpf_mul(values.b_next, values.b_next, tmp);
			// tmp = a_prev + b_prev
			mpf_add(tmp, values.a_prev, values.b_prev);
			// b_next = b_next / tmp
			mpf_div(values.b_next, values.b_next, tmp);

			// CALCULA P_NEXT
			// p_prev = p_next
			mpf_set(values.p_prev, values.p_next);
			// p_next = a_next + 1
			mpf_add_ui(values.p_next, values.a_next, 1);
			// p_next = p_next * p_prev
			mpf_mul(values.p_next, values.p_next, values.p_prev);
			// p_next = p_next * b_next
			mpf_mul(values.p_next, values.p_next, values.b_next);			
			// tmp = b_next + 1
			mpf_add_ui(tmp, values.b_next, 1);
			// p_next = p_next + tmp
			mpf_div(values.p_next, values.p_next, tmp);

		}

	// INICIO DO PROCESSO EM PARALELO
	} else {

		// Inicializa as threads
		pthread_t a_next_t;
		pthread_t b_next_t;
		pthread_t p_next_t;

		// Inicializa a interação
		for(i=0; i<IT; ++i) {
			
			// a_prev = a_next
			mpf_set(values.a_prev, values.a_next);
			// b_prev = b_next
			mpf_set(values.b_prev, values.b_next);

			
			// Executa as threads para calcular 'a' e 'b' e uma parte de 'p'
			pthread_create(&a_next_t, NULL, calc_a, &values);
			pthread_create(&b_next_t, NULL, calc_b, &values);
			pthread_create(&p_next_t, NULL, calc_p, &values);

			// espera todas as threads terminarem
			pthread_join(a_next_t, NULL);
			pthread_join(b_next_t, NULL);
			pthread_join(p_next_t, NULL);

			// Calcula o resto de 'p'. Como 'p' depende de a_next e b_next
			// não tem como executar em paralelo

			// p_next = a_next + 1
			mpf_add_ui(values.p_next, values.a_next, 1);
			// p_next = p_next * p_prev
			mpf_mul(values.p_next, values.p_next, values.p_prev);
			// p_next = p_next * b_next
			mpf_mul(values.p_next, values.p_next, values.b_next);			
			// tmp = b_next + 1
			mpf_add_ui(tmp, values.b_next, 1);
			// p_next = p_next + tmp
			mpf_div(values.p_next, values.p_next, tmp);
		}

			
	}

	// Imprime o valor de PI e o numero de casas de precisão
	int n = 50;
	gmp_printf("\n\n   PI = %.*Ff\n\n", n, values.p_next);

	return 0;
}


// Calcula a_next
void *calc_a (void *param){

	data *val;
	val = (data *) param;
	
	// Variavel para realizar calculo
	mpf_t tmp;
	mpf_init(tmp);
	
	// tmp = sqrt(a_prev)
	mpf_sqrt(tmp, val->a_prev);
	// a_next = 1 / tmp
	mpf_ui_div(val->a_next, 1, tmp);
	// a_next = a_next + tmp
	mpf_add(val->a_next, val->a_next, tmp);
	// a_next = a_next / 2
	mpf_div_ui(val->a_next, val->a_next, 2);

	pthread_exit(0);
}

// Calcula b_next
void *calc_b (void *param){
	
	data *val;
	val = (data *) param;

	// Variavel para realizar calculo
	mpf_t tmp;
	mpf_init(tmp);

	// tmp = sqrt(a_prev)
	mpf_sqrt(tmp, val->a_prev);

	// b_next = b_prev + 1
	mpf_add_ui(val->b_next, val->b_prev, 1);
	// b_next = b_next * tmp
	mpf_mul(val->b_next, val->b_next, tmp);
	// tmp = a_prev + b_prev
	mpf_add(tmp, val->a_prev, val->b_prev);
	// b_next = b_next / tmp
	mpf_div(val->b_next, val->b_next, tmp);

	pthread_exit(0);
}

void *calc_p (void *param){

	data *val;
	val = (data *) param;

	// Variavel para realizar calculo
	mpf_t tmp;
	mpf_init(tmp);
	
	// p_prev = p_next
	mpf_set(val->p_prev, val->p_next);
	
	pthread_exit(0);
}

